Resonances of dynamical checkerboard states in Josephson arrays with self-inductance 
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We study the dynamics of fully frustrated, underdamped Josephson arrays. Experiments reveal 
remarkable similarities among the dc current- voltage characteristics of several kinds of square and 
triangular arrays, where two resonant voltages are observed. Simulations indicate that a dynamical 
checkerboard solution underlies these similarities. By assuming such a solution, we reduce the 
governing equations to three coupled pendulum equations, and thereby calculate the voltages of the 
intrinsic resonances analytically. 

PACS Numbers: 74.50.+r, 05.45. +b, 74.40. +k 



Discrete arrays of nonlinear oscillators can exhibit di- 
verse spatiotemporal patterns. Examples include kinks in 
chains of coupled pendula, neuromuscular waves in the 
intestine, and modulated waves and chaos in networks 
of phase-locked loops |Q. Although such oscillator ar- 
rays are difficult to analyze completely, one can often use 
the symmetries of the system to construct simple pat- 
terns composed of spatially repeated "unit cells" . Then 
the governing equations reduce to a much smaller set of 
equations for each unit cell. This strategy has been used 
recently to construct rotating spiral waves in a model of 
discrete excitable media || . 

In this paper we present experiments, simulations, and 
analysis on a broad class of discrete arrays of nonlinear 
Josephson-junction oscillators. Networks ranging from 
single square and triangular plaqucttes to one- and two- 
dimensional arrays are considered. In each case, the junc- 
tions are identical and underdamped, and the arrays are 
driven by a dc bias current. These arrays are known to 
have complicated dynamics, and only a few analytical 
results have been obtained However, one class of 

arrays is relatively tractable: fully frustrated arrays, i.e., 
arrays subjected to an applied magnetic field of / = 1/2 
flux quanta per plaquette on average. Then the ground 
state is highly symmetric - the junction phases adopt a 
checkerboard pattern [Q. This pattern can be robust: 
even if the array is driven by an applied dc current Q , or 
a combined ac+dc current, as in studies of giant Shapiro 
steps HH, the checkerboard can persist, though now in 
the form of a propagating state. 

All of the fully frustrated arrays we have tested ex- 
perimentally exhibit strikingly similar dc current-voltage 
(/V") characteristics, with two resonant voltages V+ ~ 
^oiLsC)- 1 / 2 and V_ ~ ^(LjC)- 1 ' 2 . (Here $ is the 
flux quantum, L s is the loop inductance, C is the junc- 
tion capacitance, and Lj = $o/2irI c is the Josephson 
inductance for a junction with a critical current I c ). We 



show below that these similarities in the IV curves can 
be explained by assuming that the junction phases are 
organized into a dynamical checkerboard solution. This 
ansatz reduces the governing circuit equations to three 
coupled pendulum equations, which in turn allows us to 
obtain analytical predictions for the observed resonant 
voltages. 

Previous authors have used the ansatz of a dynami- 
cal checkerboard, or some other symmetric pattern, to 
obtain reduced equations for frustrated Josephson ar- 
rays Our analysis extends this work in three 
respects: (1) Inductance effects are essential to explain 
some of the experimental results. Therefore, we include 
self-inductance in the analysis and simulations as a first 
approximation. (2) The governing equations for the dy- 
namical checkerboard are shown to have the same alge- 
braic form for all eight arrays in Figure ^ up to topolog- 
ical factors that encode the array geometry. This shared 
structure allows a unified analysis of the different net- 
works. (3) Intrinsic resonant voltages are predicted ana- 
lytically from the reduced equations. 

Figure [l] depicts the arrays and their dynamical 
checkerboard states. Junctions (not shown) exist along 
each branch, and current is applied at each node along 
the top and bottom edges. Throughout this paper, "ver- 
tical" refers to the direction of current injection. We 
consider junctions with identical critical currents I c , re- 
sistances R n and capacitances C. Experiments have been 
performed on arrays of Nb-A^O^-Nb junctions Jl(| with 
small damping T — f3 c ' , where 5 < (3 C < 25 is the Mc- 
Cumber parameter. The magnetic field is applied nor- 
mal to the arrays by an external solenoid. The two- 
dimensional penetration depth A = Lj/L s , which can 
be viewed as a measure of the discreteness of the ar- 
rays, ranges from 0.5 < A < 2. By use of a diagnos- 
tic technique pi], we determine the parameters for each 



1 



array: first, the measurement of R n and I C (T), with 
I c (0)R n = 1.9 mV, yields the temperature-dependent 
Josephson inductance Lj(T); second, from the Fiske 
modes of the diagnostic we obtain C, L s , and the mu- 
tual inductance for nearest neighbors, which we find to 
be small |l^Jl^]. To facilitate comparisons between dif- 
ferent types of arrays, the triangular arrays were built by 
decimating every other horizontal junction from the cor- 
responding square array; thus their self-inductances are 
equal. 

Figure ^| shows the measured IV curves (I is the cur- 
rent per vertical junction normalized by I c , and V is the 
voltage per row) for three different geometries at / = 1/2. 
The signature of all these IVs is the appearance of jumps 
at two resonant voltages, V + and V-. The upper step, 
which ends at V+, is nearly vertical and independent of 
temperature. For 1HSQ and 2DSQ arrays (the arrays we 
have studied most extensively) , V+ is also independent of 
the number of cells in the x-direction . These results 
suggest that local geometrical properties, i.e. L s , de- 
termine the voltage. A temperature-independent voltage 
V+ ~ $0/ VL S C is expected from dimensional arguments 
and is found to be approximately correct experimentally 
and in simulations. When varying the magnetic frustra- 
tion, V+ has the usual periodicity in / and reflection sym- 
metry about / = 1/2 but it is found to be almost inde- 
pendent of / for 0.2 < / < 0.5. Thus, although the value 
of V+ is similar to the Eck peak in ID parallel arrays Jli) , 
it does not have the sine-like /-dependence observed in 
that case |0|. Moreover, this upper step does not ap- 
pear for / < 0.2 or for A < 0.5 as the system switches 
directly into row-switched states. The lower voltage V- , 
on the other hand, is temperature-dependent, suggest- 
ing a dependence on the Josephson inductance; namely, 
V- ~ $0/ \/Lj(T)C, up to a factor of order unity. 

The dynamical origin of these two resonances is re- 
vealed through numerical simulation of the complete ar- 
rays. The governing equations result from current conser- 
vation at each node and from fluxoid quantization around 
each loop. The current Ij through each junction, in the 
resistively and capacitively shunted junction model, is 
Ij = sin (f>j + Tcj)j + ipj where the currents are in units of 
I c , and time is in units of (L/C) 1 / 2 , the inverse of the 
plasma frequency. Fluxoid quantization demands that 
<j)j = —2nf — icir/A, where the sum is around a loop, 
7 c i r is the loop current, and only self-(loop)-inductances 
are considered. Simulated IVs for square arrays are con- 
sistent with the experiments [[l5| . 

Furthermore, the simulations at / = 1/2 suggest that 
solutions with wavelength equal to two plaquettes, as 
shown in Figure^, underlie the observed numerics. These 
dynamical checkerboard solutions are constituted by sin- 
gle plaquettes (SQ and TR as shown) paired with their 
symmetric counterparts. Under these symmetry con- 
straints, the number of relevant variables is greatly re- 
duced and the dynamics of the whole array is governed 



by a set of three coupled equations. Moreover, the equa- 
tions for each of the eight arrays in Figure [j] can be recast 
in the same algebraic form, up to constant factors that 
encode the array topology. In this unified formulation, 
the fluxoid quantization condition at / = 1/2 is 

4>i - 02 - h(j>3 = -it + h/(\m y ) (1) 
and the two current conservation conditions arc 



I\ + I2 — 21 and I2 — I\ = 2m x I^ jm 



(2) 



Here h is the number of horizontal junctions per cell (h = 
2 for square and h = 1 for triangular arrays), and m y — 
1 + Vy/2 where v y is the number of neighboring cells in 
the ^/-direction; likewise, m x = 1 + v x /2. For example, 
the 1HSQ array has h — 2, m y = 1, and m x = 2, and the 
1VTR array has h = 1, m y = 2, and m x = 1. 

It is convenient to introduce three new variables: 
= (4>i + (fa)/ 2, the average of the vertical phases; 
9 = (4>2 — 0i — 7r )/2, which measures how much the dif- 
ference of the vertical phases differs from 7r; and ip = 03, 
the top horizontal phase. When the arrays have horizon- 
tal junctions, i.e. h 7^ 0, Eqs. (0)— (@) can be rewritten as 
a system of three coupled nonlinear pendulum equations: 

(j> + T<j) - sin 9 sin = 7 (3) 
6 + T6 + cos cos 9 = -Xm x (29 + hip) (4) 
<p + r^i + simp = —Xm y (29 + hip). (5) 

Numerical simulations indicate that, for the parame- 
ter regimes considered, is approximately a uniformly 
whirling phase with frequency u>, whereas both 9 and ip 
librate sinusoidally at to. Thus, we approximate the solu- 
tion as: = cut + k, 9 = B sin (urt + a), and ip — A sintot, 
where a, k, A, B and lo are constants to be determined 
for each set of driving conditions {/, T, A}. A harmonic 
balance calculation with this assumed solution yields an 
IV characteristic with two resonant peaks, as shown by 
the dashed line in Figure I for the 1HSQ device. These 
peaks are typically associated with a loss of stability of 
the underlying dynamical state, thus explaining the ex- 
perimental steps. This semi-analytic TV approximates 
the voltages of the resonant steps reasonably well, al- 
though it does not provide a quantitative fit to the mea- 
sured currents. The details of this calculation will be 
given elsewhere |Q. 

Approximate formulas for the resonant voltages can 
be obtained via a further simplification. Since the libra- 
tion amplitudes A and B are observed to be small for a 
wide range of bias current, we regard both ip and 9 as 
small oscillations driven by the whirling mode 0. Then 
Eqs. (@)-@ can be linearized and, in the limit of small 
damping (V <C 1), these two equations have two reso- 
nance frequencies at 



u]_ = 1/2 + \m x a ± (1/2 + Am^rr) — 2Xm x 



1/2 



(6) 
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where a = 1 + hm y j (2m x ). These frequencies agree with 
the location of the peaks of the IV from the full har- 
monic balance calculation. To clarify the physical mean- 
ing of the two resonances u>± , consider the limit of small 
inductance, A > 1. Then, lo + ^/2\m x a = b + 
and uj- ~ cr -1 / 2 = b-. The corresponding voltages are 
V+ » b+ $ /(2^VTO) and U_ » 6_ * /(27r A /L J (T)C), 
as guessed earlier from dimensional arguments. Thus, 
V+ is temperature-independent while V_ depends on T; 
and 6+ and b- are combinations of topological factors of 
order unity. (For the single triangle (TR), b + — y/3, as 
suggested by Yukon and Lin ||). 

Figure || shows good agreement between the temper- 
ature dependence of the experimental (squares and tri- 
angles) and predicted (lines) normalized voltages V + /b + 
and for 1HSQ, 2DSQ and 1HTR arrays. The fact 

that Eq. (||) (dashed line) consistently overestimates the 
value for V+/b+ can be attributed to an underestimation 
of the inductance since simulations with self-inductance 
alone (pluses in the figure) agree well with the theory. 
When simulations with the full inductance matrix are 
performed (crosses), the resonant voltage moves down 
toward the experimental data. The main effect of the 
mutual inductive coupling can be approximated by an 
effective inductance L S:C g — L s {1 + M(v x + v y )} where 
M ps 0.12 is the ratio of the magnitude of the nearest- 
neighbor inductance to the self-inductance in square ar- 
rays p^ ]. The solid line corresponds to V+/b+ calculated 
from Eq. (^) with this effective inductance and predicts 
well both the full inductance simulations and the exper- 
imental data. In addition, the values and temperature 
dependence of U_/6_ agree well with the theory and 
simulations and, as expected, no inductance effects are 
visible. 

Both the horizontal junctions and the inductance play 
essential roles in the phenomena described above; with- 
out them, one or both of the observed resonances would 
be lost. To see this, consider three well-studied limit- 
ing cases. First, if the horizontal junctions are absent, 
the SQ and TR arrays are simply two-junction SQUIDS 
with inductances, and 1HSQ and 1HTR become ID par- 
allel arrays, which can be regarded as discrete versions 
of long Josephson junctions. The dynamics of these sys- 
tems are then governed by only two equations, Eqs. (|^) 
and (0) with h = 0, since p is no longer a valid dynamical 
variable. This system has a single linear resonance, at a 
frequency lu = \/2Xm x = Lu+\h=o- Second, in the limit 
A — ► oo where inductances are neglected, there is no di- 
vergence in the system since 9 = —hip/2. Thus, 
Eqs. (|]) and rtq) are combined to eliminate A (29 + hp) 
and the dynamics of the system is then governed by only 
two equations: Eq. (|^) and a reduced equation for 9, 

e + Ye+--Lm(\e] +^ cob* cob 4 =0. (7) 
cr 2 I \h ) m x J 

(These equations include as a special case the 2DSQ sys- 



tem previously studied in ^|JJ] . Note also that when A — > 
oo, the SQ and 2DSQ equations are identical but this is 
no longer true when inductances are taken into account). 
In conclusion, in the absence of inductance there is a 
single linear resonant frequency lu = a^ 1 !" 1 = u}—\x—kx>- 
Third, if there are neither horizontal junctions nor in- 
ductances, the sole dynamical equation (^|) has no reso- 
nances. 

These limiting cases also suggest a heuristic explana- 
tion for the origin of the resonances. At V_, as the 
checkerboard slides across the array, it strongly excites 
spin waves, or equivalently, "ringing" oscillations of the 
junctions at an eigenfrequency close to their character- 
istic plasma frequency (Lj(T)C)~ 1 ^ 2 . This resonance 
requires horizontal junctions but not inductance, and is 
temperature-dependent. In contrast, the resonance at U+ 
is due to excitation of electromagnetic modes of the array. 
This resonance requires inductance but not horizontal 
junctions, and depends only on local geometrical proper- 
ties of the array that are temperature- independent. The 
corresponding eigenmode is mainly related to oscillations 
in the induced flux per cell. 

The role of the horizontal junctions can be further 
explored by considering an anisotropic network where 
horizontal and vertical junctions are fabricated with the 
same process and differ only in their area. Then the ra- 
tio of their critical currents, r\ = J c i/J C 3, parametrizes 
the anisotropy of the array. In this case, the results 
from equations are still valid with a renormal- 

ized m'y — m y 77, which specifies the participation of the 
horizontal junctions in the dynamics of the array. 

Finally, we emphasize that our analysis is based on the 
assumption of a dynamical checkerboard state. Unfortu- 
nately, very little is known about the conditions for its 
global stability. If it is unstable, or if the array orga- 
nizes itself in some alternative stable state, the dynam- 
ics are not yet understood. For instance, simulations of 
2DTR arrays with low j3 c seem to show a different solu- 
tion (the "ribbon state") where the horizontal junctions 
are essentially inactive and a checkerboard is formed by 
double cells J^] . This solution is analogous to the striped 
columnar dynamical configuration observed in / = 1/2 
Shapiro steps in square arrays pq| . In both cases, the 
arrays effectively behave as a collection of in-phase rows 
and, thus, have only one resonant voltage at V+|h=o- Ex- 
periments on triangular arrays will address these issues 
separately p7| . Also, if the junctions are highly under- 
damped, the checkerboard state in 2DSQ arrays can slide 
chaotically B. The conditions for the stability and the 
temporal periodicity of the checkerboard state, and the 
dynamics associated with other possible states, are chal- 
lenging problems for future investigation. 
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FIG. 1. Arrays and their checkerboard states: Single 
square (SQ) and single triangular (TR) plaquettes; ID square 
(1HSQ and 1VSQ) and ID triangular (1HTR and 1VTR) 
arrays; 2D square (2DSQ) and 2D triangular (2DTR) ar- 
rays. The phases of the vertical junctions <f>i and 2 and 
the horizontal junction </>3 are defined in the SQ and TR di- 
agrams. The checkerboard solutions are constructed by com- 
bining these single plaquettes with their symmetric counter- 
parts. 
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FIG. 2. Experimental IV curves for three arrays: 1HTR 
array (1x9 plaquettes) with /3 C = 8 and A = 0.64; 1HSQ 
array (1 x 7) with f3 c = 11 and A = 0.76; and 2DSQ array 
(7 x 7) with (3 C = 20 and A = 0.92. Dashed line, IV from 
harmonic balance for the 1HSQ array with the same f3 c and 
an effective A e fi = 0.61 which accounts for mutual inductance 
effects. V+ and V- are indicated for 1HTR. 
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FIG. 3. Normalized resonant voltages vs. T, measured 
for the arrays presented in Fig. |^: 1HTR (solid triangles), 
1HSQ (solid squares), and 2DSQ (open squares). Dashed 
lines, Eq. (^|). Solid lines, Eq. ([]) with L s , e fi as defined in 
the text. Values of V+/6+ from numerical simulations of the 
1HSQ array with only self-inductance (pluses) and the full 
inductance matrix (crosses) are also shown. The 1HTR array 
was built from a 1HSQ by decimation of every other horizontal 
junction. Therefore, the geometrical unit cells are identical 
and L s — 23 pH and C = 300 fF for the three arrays. 
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